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ABSTRACT 


We present JWST Early Release Science (ERS) coronagraphic observations of the super-Jupiter 
exoplanet, HIP 65426 b, with the Near-Infrared Camera (NIRCam) from 2—5 yum, and with the Mid- 
Infrared Instrument (MIRI) from 11-16 ym. At a separation of ~0.82” (87739° au), HIP 65426 b is 
clearly detected in all seven of our observational filters, representing the first images of an exoplanet 
to be obtained by JWST, and the first ever direct detection of an exoplanet beyond 5 wm. These 
observations demonstrate that JWST is exceeding its nominal predicted performance by up to a factor 
of 10, with measured 5c contrast limits of ~4x10~° (~2.4 wJy) and ~2x1074 (~10 pJy) at 1” for 
NIRCam at 3.6 wm and MIRI at 11.3 um, respectively. These contrast limits provide sensitivity to 
sub-Jupiter companions with masses as low as 0.3Mjyp beyond separations of ~100 au. Together 
with existing ground-based near-infrared data, the JWST photometry are well fit by a BT-SETTL 
atmospheric model from 1—16 wm, and span ~97% of HIP 65426 b’s luminous range. Independent of 
the choice of forward model atmosphere we measure an empirical bolometric luminosity that is tightly 
constrained between log(Lpo1/Lo)=—4.35 to —4.21, which in turn provides a robust mass constraint 
of 7.1£1.1 Mjyp. In totality, these observations confirm that JWST presents a powerful and exciting 
opportunity to characterise the population of exoplanets amenable to direct imaging in greater detail. 


1. INTRODUCTION 


Across the last twenty-five years a variety of obser- 
vational techniques have been employed to uncover and 
characterise the current population of over 5000 con- 
firmed exoplanets (NASA Exoplanet Science Institute 
2020). Of these techniques, the direct detection of pho- 
tons from an exoplanetary atmosphere — direct imaging 
— remains one of the most challenging due to the sub- 
stantial contrast in flux between host stars and their exo- 
planetary companions. The emitted flux of an exoplanet 
can be many magnitudes fainter than the stellar diffrac- 
tion halo at its angular separation, and bespoke instru- 
mentation (e.g., Macintosh et al. 2014; Beuzit et al. 
2019) and/or image post-processing (e.g., Chauvin et al. 
2004; Marois et al. 2008) is needed to isolate the exo- 
planet emission. Nevertheless, these “high contrast” ob- 
servations offer considerable advantages to other tech- 
niques. At present, direct imaging is the most viable 
technique for characterising the population of exoplan- 
ets at orbital separations greater than ~10 au. Further- 
more, beyond the large population of irradiated, close- 
in planets with atmospheric measurements obtained via 
exoplanet transit observations, high contrast observa- 
tions of exoplanet emission remain the most readily ac- 
cessible path towards the characterisation of exoplanet 
atmospheres. Constraints on atmospheric composition 
may improve our understanding of exoplanet formation 
and evolution (e.g., Oberg et al. 2011; Morley et al. 
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2019; Zhang et al. 2021; Molliére et al. 2022), although 
these determinations can be highly dependent on the 
post-formation accretion of solid material. Compared 
to close-in transiting exoplanets, directly imaged plan- 
ets present a distinct advantage in this regard, as they 
are easier to detect at younger ages where they are less 
likely to have experienced significant migration and/or 
accretion. Additionally, at young ages bulk properties 
such as temperature, radius, and bolometric luminosity 
provide independent constraints on formation conditions 
(Marley et al. 2007; Marleau & Cumming 2014) that 
can be contrasted to atmosphere-driven conclusions on 
formation. Finally, the study of exoplanet atmospheres 
continues to advance towards smaller, and more Earth- 
like exoplanets, and could ultimately lead to the discov- 
ery of life outside our Solar System (Schwieterman et al. 
2018). 

At present ~20 planetary mass companions (PMCs) 
have been detected and characterised through high con- 
trast observations (Currie et al. 2022), and all exoplanets 
directly imaged to date have estimated or dynamically- 
measured masses 22 Mjy). Despite the small sample 
size, this subset of objects has driven significant de- 
velopments in our overall understanding of exoplanet 
atmospheres and architectures. For example, the for- 
mation of exoplanets through gravitational instability is 
rare (Vigan et al. 2017; Nielsen et al. 2019; Vigan et al. 
2021); brown dwarfs and exoplanets may exhibit differ- 
ent eccentricity distributions (Bowler et al. 2020); clouds 
and disequilibrium chemistry influence measured spec- 
tra (e.g., Skemer et al. 2012; Apai et al. 2013; Morley 
et al. 2014); similar to brown dwarfs, exoplanets can be 
variable (Zhou et al. 2019, 2020; Vos et al. 2022); wide 
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separation exoplanets may be more prevalent within 
systems that also host a circumstellar disk (Meshkat 
et al. 2017); and transits of exomoons may be observable 
around isolated planetary mass objects (Limbach et al. 
2021). Nevertheless, upgraded or new observatories and 
instruments (e.g., Gardner et al. 2006; Males et al. 2018; 
Chilcote et al. 2020) will be necessary to directly image 
and characterise a wider diversity of closer separation 
and/or lower mass exoplanets at a higher precision and 
across a broader wavelength range. 


1.1. High Contrast Observations with JWST 


Launched on December 25th, 2021, JWST (Gard- 
ner et al. 2006) is an international collaboration be- 
tween NASA, ESA, and the CSA, and the first large 
strategic mission of the NASA Astrophysics Division 
to launch since the 1990’s (National Academies of Sci- 
ences & Medicine 2017). With a primary mirror di- 
ameter of 6.5 m, an operational wavelength range 
from 0.6—28.1 wm, and a diverse range of instrumental 
modes, JWST presents a revolutionary opportunity for 
scientific exploration and discovery across all branches 
of astronomy and astrophysics. Within this remit, high 
contrast observations of exoplanets and exoplanetary 
systems are no exception. 

JWST is located at the second Sun—Earth Lagrange 
point, far from the thermal background, telluric contam- 
ination, and wavefront aberrations generated by Earth’s 
atmosphere. The combination of excellent optical per- 
formance (~75—130 nm RMS wavefront error, depend- 
ing on the instrument), highly stable wavefront (<2 nm 
drift over a few hours), and large telescope aperture, 
enables JWST to reach better photometric and spec- 
troscopic limiting sensitivities than both past or current 
ground- (i.e., 8—10 m class telescopes) and space-based 
(e.g., Hubble, Spitzer) observatories (Rigby et al. 2022). 
It is not only this increased sensitivity that improves our 
ability to detect and characterise faint objects such as 
exoplanets, but its combination with JWST’s access to 
the near- and mid-infrared. At these wavelengths, the 
flux emitted from a hotter host star steadily decreases 
as a function of increasing wavelength, whereas the flux 
emitted from cooler exoplanetary companions reaches a 
peak. Hence, the natural star-planet contrast is dimin- 
ished. To realise these advantages, JWST offers a selec- 
tion of instrumental modes designed for, or that can be 
applied to, high contrast observations. Specifically, both 
NIRCam (Rieke et al. 2005) and MIRI (Rieke et al. 2015) 
are equipped with coronagraphic masks (Krist et al. 
2009; Boccaletti et al. 2015), NIRISS (Doyon et al. 2012) 
is equipped with a non-redundant mask which enables 
aperture masking interferometry (AMI; Sivaramakrish- 


nan et al. 2012), and although lacking any hardware for 
starlight suppression, both NIRSpec (Bagnasco et al. 
2007) and MIRI are equipped with integral field units 
(IFUs; Wells et al. 2015; Boker et al. 2022). 

In anticipation of the range of capabilities that JWST 
would provide, a similar range of predictions and simula- 
tions were constructed in an effort to forecast its poten- 
tial for exoplanet imaging science. Brande et al. (2020) 
predict that Saturn and Jupiter mass planets should be 
detectable with MIRI coronagraphy at 1—5 au sepa- 
rations across a sample of nearby (S10 pc) M-dwarfs. 
Similarly, Carter et al. (2021a) demonstrate that both 
NIRCam and MIRI coronagraphy may be sensitive to 
sub-Jupiter mass planets for the majority of stars within 
the 2 Pictoris moving group (GPMG) and TW Hydrae 
Association (TWA), albeit at separations >20—50 au. 
Ray et al. (submitted) further expands on the work 
of Carter et al. (2021a), and shows that sub-Jupiter 
mass planets could be detectable from 1—20 au for sev- 
eral stars within GPMG, TWA, and the Taurus Auriga 
Association with NIRISS AMI. At even closer separa- 
tions, Beichman et al. (2020) predict that MIRI coro- 
nagraphy of a Centauri A may be sensitive to ~5 Re 
companions from 0.5—2.5 au. Danielski et al. (2018) 
demonstrate that already discovered companions (e.g., 
HR 8799b/c/d Marois et al. 2008, 2010), 8 Pictorisb 
(Lagrange et al. 2010)) should also be detectable with 
MIRI coronagraphy, and NH3 absorption could be iden- 
tified for a subset of targets. Finally, Patapis et al. 
(2022) show that molecular mapping is possible with 
the MIRI IFU and may result in the detection of atmo- 
spheric species such as HzO, CO, CO2, CH4, NH3, and 
PH3. 

These predictions are based on ground-based testing 
and observatory simulations, however, the true capabil- 
ities of JWST hinge on its on-sky performance. Prelim- 
inarily, the performance of both the NIRCam and MIRI 
coronagraphic modes exceeded expectations during ob- 
servatory commissioning (Girard et al. 2022; Kammerer 
et al. 2022; Boccaletti et al. 2022), but the first scien- 
tific demonstrations of JWST’s capabilities are being 
conducted as part of the Director’s Discretionary Early 
Release Science (ERS) Programs!. Our ERS program 
“High Contrast Imaging of Exoplanets and Exoplanetary 
Systems with JWST” (ERS-01386, Hinkley et al. 2022a) 
is the only ERS program that has tested the high con- 
trast exoplanet imaging modes of JWST and includes: 
coronagraphic imaging from 2—16 jm of the known ex- 
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oplanet HIP 65426b (Chauvin et al. 2017, this work) 
and circumstellar disk HD 141569 A (Weinberger et al. 
1999, Millar-Blanchaer et al. in preparation, Choquet 
et al. in preparation), spectroscopy from 1—28 ym of 
the PMC VHS J125601.92-125723.9 ABb (VHS 1256b, 
Gauza et al. 2015, Miles et al. submitted), and AMI ob- 
servations of HIP 65426 at 3.8 wm (Ray et al. in prepara- 
tion). This program is rapidly disseminating these cru- 
cial initial data, and demonstrating the true capabilities 
of JWST for high contrast imaging and spectroscopy for 
the first time. Furthermore, we will provide a range of 
science enabling products (e.g., data analysis pipelines, 
recommendations for best observing practices) to the 
community to support their own proposals and investi- 
gations in Cycle 2 and beyond (Hinkley et al. 2022a). 

In this work we focus exclusively on the coronagraphic 
imaging observations of the HIP 65426 system within 
this ERS program, and their context within a broader 
understanding of JWST as a tool for high contrast imag- 
ing. 


1.2. HIP 65426b 


Discovered by Chauvin et al. (2017), HIP 65426 b is a 
super-Jupiter mass exoplanet at a wide physical sep- 
aration of 110733 au (Cheetham et al. 2019) to the 
star HIP 65426 (A2V, 2MASS J=6.826, AJ—K=0.055, 
Mo=1.9640.04). HIP 65426 is located at a distance 
of 107.49+0.40 pe (Gaia Collaboration 2020), has no 
signs of binarity from radial velocity and sparse aperture 
masking observations (Chauvin et al. 2017; Cheetham 
et al. 2019), and is a fast rotator (vsin(i)=299+9 km st, 
Chauvin et al. 2017). Furthermore, HIP 65426 is a likely 
member of the Lower Centaurus-Crux association as de- 
rived from its proper motion and radial velocity mea- 
surements (89% probability, Gagné et al. 2018), con- 
straining its age to 1444 Myr. This association has 
grown in interest over time with the increasing number 
of directly imaged exoplanet discoveries within it (e.g., 
HD 95086 b Rameau et al. 2013, PDS 70b/c Keppler 
et al. 2018, and TYC 8998b/c Bohn et al. 2020. 

Although HIP 65426 b was initially observed with a 
combination of low resolution spectroscopy and pho- 
tometry from ~1—2 wm (Chauvin et al. 2017), follow 
up observations expanded this coverage to ~1—5 pm 
(Cheetham et al. 2019; Stolker et al. 2020b), includ- 
ing a medium-resolution spectrum from ~2—2.5 ym 
(R~5500, Petrus et al. 2021). Photometric analy- 
sis has demonstrated that HIP 65426b is similarly lo- 
cated to mid-to-late L-dwarfs in colour-magnitude di- 
agrams (Figure 1), and lies between already discov- 
ered early L-dwarf exoplanet companions (e.g., 3 Picb, 
HD 106906b) and those at the L/T transition (e.g., 


M, [mag] 
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Figure 1. Colour-magnitude diagram showing the posi- 
tion of HIP 65426 b (Chauvin et al. 2017) relative to both 
the population of brown dwarf objects (circles) and a selec- 
tion of directly imaged planetary mass companions (PMCs, 
hexagons), as obtained from Best et al. (2020). 


HR 8799 c/d/e). Using combined photometric and spec- 
troscopic observations, Petrus et al. (2021) performed an 
atmospheric forward model analysis of HIP 65426 b, in- 
dicating that it has T.g¢=1560+100 K, log(g)<4.40 dex, 
[M/H]=0.0510'35 dex, and the atmospheric carbon-to- 
oxygen ratio has an upper limit of C/O<0.55.  Fur- 
thermore, Petrus et al. (2021) also detect the presence 
of H2O and CO in the atmosphere of HIP 65426b us- 
ing cross-correlation molecular mapping, in addition to 
non-detections of CH, and NH3. Finally, from evolu- 
tionary model analyses to the data reported in Chauvin 
et al. (2017), Marleau et al. (2019) estimate the mass 
of HIP 65426 b to be 9.9t}'3 Mjup or 10.943°5 Myup for 
hot- or cold-start initial entropy conditions, respectively 
(see e.g., Marley et al. 2007; Spiegel & Burrows 2012). 


In Section 2 we describe our JWST observations of 
HIP 65426b and all necessary data reduction. In Sec- 
tion 3 we describe the analysis steps taken to produce 
residual starlight subtracted images and measurements 
of contrast performance. We present a discussion of 
these observations in the context of both the overall per- 
formance of JWST and our individual understanding of 
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HIP 65426b in Section 4. Finally, we summarise our 
conclusions in Section 5. 


2. OBSERVATIONS & DATA REDUCTION 


The NIRCam and MIRI coronagraphic imaging ob- 
servations of HIP 65426b presented here were taken as 
part of program ERS-01386 (Hinkley et al. 2022a) and 
exist as a subset of a broad range of observations to as- 
sess the performance of JWST’s high contrast imaging 
and spectroscopic modes with respect to the study of 
exoplanetary systems (Hinkley et al. 2022a). 


2.1. Observational Structure 


The observational strategies used for this program 
were adopted following the recommended best practices 
as known prior to launch and described in the JWST 
user documentation?’?. All observations of HIP 65426 b 
are repeated at two independent roll angles separated 
by ~10° to enable subtraction of the residual stellar 
point spread function (PSF) through angular differen- 
tial imaging (ADI, Miiller & Weigelt 1985; Liu 2004; 
Marois et al. 2006). Although a large number of rolls 
across a larger angular range would be desirable for an 
optimal subtraction using this technique, the combina- 
tion of lengthy exposure times, increased overheads, and 
spacecraft orientation constraints prohibit an observing 
strategy more complex than described. Given the max- 
imum possible roll offset for JWST at any given epoch 
is 14°, a larger roll offset than we have adopted would 
also require multi-epoch observations. 

We also perform similar observations of a bright 
reference star, HIP 68245 (B2IV, 2MASS J=4.628, 
AJ—K=0.137), to additionally enable subtraction 
through reference differential imaging (RDI). This star 
was selected as it is: a) bright, therefore reducing the 
exposure time required to attain a similar signal-to-noise 
as the science target, b) a similar spectral type to the 
science target, therefore reducing the impact of spectral 
type mismatch, c) is relatively close (~10°) to the sci- 
ence target, therefore reducing slew overheads and min- 
imising position dependent wavefront drift between sci- 
ence and reference observations, and d) has no evidence 
of binarity as determined by VET /SPHERE AMI ob- 
servations (Proposal ID: 108.22CD). For information on 
selecting a suitable reference star, see the JWST User 


2 https: //jwst-docs.stsci.edu/jwst-near-infrared-camera/nircam- 
observing-strategies /nircam-coronagraphic-imaging- 
recommended-strategies 

3 https: //jwst-docs.stsci.edu/jwst-mid-infrared-instrument /miri- 
observing-strategies /miri-coronagraphic-recommended- 
strategies 


Documentation*. All exposure settings for the science 
and reference observations are shown in Table 1. 

An RDI based subtraction is likely to reach superior 
contrast limits to ADI from pre-launch predictions (La- 
joie et al. 2016; Perrin et al. 2018; Carter et al. 2021b), 
and is therefore also more representative of the optimal 
performance of JWST coronagraphy (also see Section 
4.1). The exposure settings for these reference observa- 
tions were chosen to reach an approximately equivalent 
fraction of full well saturation per integration as the cor- 
responding target observations. Additionally, each ref- 
erence observation was repeated at nine separate dither 
positions following small-grid dither patterns 9-POINT- 
CIRCLE and 9-POINT-SMALL-GRID for the NIRCam 
and MIRI observations, respectively (Soummer et al. 
2014; Lajoie et al. 2016). The goal of this strategy is to 
produce a small library of reference PSFs for each sci- 
ence exposure which captures different misalignments 
between the star and the center of the coronagraphic 
mask, and can in turn facilitate more advanced PSF sub- 
traction techniques (e.g., KLIP, Soummer et al. 2012). 
Further discussion on the relative benefits between ADI 
and RDI subtraction strategies, or a combination of the 
two, with respect to these JWST observations can be 
found in Section 4.1. 

For MIRI, we also add background observations to 
both our science and reference observations in both fil- 
ters of a nearby “empty” region of sky (as identified 
in WISE images, Wright et al. 2010) separated ~1.5' 
away the target star to measure the stray light “glow 
stick” that is inherent to MIRI coronagraphic observa- 
tions (Boccaletti et al. 2022). Specifically, this posi- 
tion corresponds to a right ascension and declination of 
@ =13°24'44.2915", 6=—51°29'31.54”, respectively. To 
best match the science and reference observations, the 
exposure parameters for each background observation 
exactly match the parameters for a single roll/dither of 
the associated science or reference target. These ob- 
servations were intended to be performed at two sep- 
arate dither positions to identify astrophysical sources 
that might impact the background subtraction, however, 
due to a previously unresolved issue they were instead 
repeated at an identical pointing (Dean Hines, private 
communication). 

The NIRCam and MIRI observations were executed 
as two separate non-interruptible sequences, ensuring 
that observations between rolls, and also between sci- 
ence, reference, and background targets, are minimally 


4 https: //jwst-docs.stsci.edu/methods-and-roadmaps/jwst- 
high-contrast-imaging/jwst-high-contrast-imaging-proposal- 
planning/hci-psf-reference-stars 
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Filter Amean (uM) Wer (zm) Mask Readout Neroups Nints  texp (8)  Naithers Nrons  ttotal (8) 
HIP 65426 
F250M 2.523 0.179 MASK335R DEEP8 15 4 1235.892 1 2 2471.784 
F300M 3.067 0.325 MASK335R DEEP8 15 4 1235.892 1 2 2471.784 
F356W 3.580 0.769 MASK335R DEEP8 15 2 617.946 1 2 1235.892 
F410M 4.084 0.436 MASK335R DEEP8 15 2 617.946 1 2 1235.892 
F444W 4.397 0.979 MASK335R DEEP8 15 2 617.946 1 2 1235.892 
F1140C 11.307 0.608 FQPM1140 FASTRI 101 41 1002.102 1 2 2004.204 
F1550C 15.514 0.703 FQPM1550 FASTRI 250 60 3609.341 1 2 7218.682 
HIP 68245 
F250M 2.523 0.179 MASK335R. MEDIUM8 4 4 166.852 9 if 1501.669 
F300M 3.067 0.325 MASK335R MEDIUM8 4 4 166.852 9 1 1501.669 
F356W 3.580 0.769 MASK335R MEDIUM8 4 2 83.426 9 1 750.835 
F410M 4.084. 0.436 MASK335R, MEDIUM8 4 2 166.852 9 1 750.835 
F444W 4.397 0.979 MASK335R MEDIUM8 4 2 83.426 9 1 750.835 
F1140C 11.307 0.608 FQPM1140 FASTRI 52 10 126.791 9 1 1141.116 
F1550C 15.514 0.703 FQPM1550 FASTRI 100 19 459.706 9 1 4137.356 


Table 1. Target and reference exposure settings. Background observations were also performed for the MIRI F1140C and 
F1550C filters with parameters identical to two exposures of a single roll or dither of the target and reference observations, 
respectively (see Section 2.1). Filter mean wavelengths (Amean) and bandwidths (We) are taken from spaceKLIP (see Section 


2.2). See https: //jwst-docs.stsci.edu/understanding-exposure-times for further detail on JWST exposure settings. 


separated in time. This reduces the extent to which the 
wavefront can vary across observations, due to variations 
in the telescope mirror alignment, the thermal evolution 
of the telescope, or both (Perrin et al. 2018; Carter et al. 
2021b). Changes in the wavefront will lead to variations 
in the residual PSF between exposures, hinder our abil- 
ity to perform an optimal subtraction of these residual 
PSFs, and suppress the overall achievable contrast. 


2.2. spaceK LIP 


For all observations we perform data reduction us- 
ing the newly developed and publicly available python 
package, spaceKLIP® (Kammerer et al. 2022). Briefly, 
spaceKLIP takes a collection of data products from the 
jwst pipeline®:’ (Bushouse et al. 2022) as inputs, and 
generates PSF subtracted images, contrast curves, and 
measurements of companion photometry and astrome- 
try. The majority of this functionality is provided by 
the underlying pyKLIP (Wang et al. 2015) package, with 
spaceKLIP providing a user friendly interface, stream- 
lined code execution, custom JWST data reduction rou- 
tines, and built-in plotting procedures. 


5 https: //github.com/kammerje/spaceKLIP 
6 https: //jwst-pipeline.readthedocs.io 


7 All data were processed using pipeline version, CAL_VAR=1.6.2, 


and calibration reference data, CRDS_CTX=jwst_0942.pmap 


2.3. NIRCam Coronagraphy 


The NIRCam observational sequence was executed 
from 23:00 July 29th to 05:16 July 30th 2022 UTC, with 
exposures taken for HIP 65426 b using the MASK335R 
round coronagraphic mask (Krist et al. 2010) in the 
F250M, F300M, F410M, F356W, and F444W filters at 
one roll angle, then once again at a second roll angle, and 
then finally for the reference star. This sequence struc- 
ture significantly reduces overheads, as once the target 
acquisition has been performed it is not necessary to 
reacquire the target to switch the observational filter. 

We begin data reduction using the Stage 0 (*un- 
cal.fits) files as generated by the jwst pipeline. These 
products are then processed to Stage 1 (*rateints.fits) 
files using spaceKLIP, which follows a slightly modified 
version of the jwst pipeline. Where possible the JWST 
detectors have reference pixels that can be used to track 
pixel drifts due to readout electronics. In the absence of 
such a reference, these variations may be misinterpreted 
as “jumps” in the up-the-ramp (MULTIACCUM) de- 
tector readout. The NIRCam coronagraphic subarrays 
do not have any embedded reference pixels and default 
pipeline processing leads to multiple erroneous jump de- 
tections and greatly increased noise in processed images. 
Therefore, we manually define all pixels within a four 
pixel border of each image as reference pixels within the 
pipeline to mitigate these effects. Additionally, we iden- 
tify a significant improvement in image quality by skip- 
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ping the dark current subtraction step that is turned on 
as a default in the jwst pipeline. This improvement is 
primarily driven by the poor quality of the dark current 
reference file, which exhibits a large number of hot pix- 
els, persistence, and cosmic rays. Once this reference file 
is improved through further calibration observations, we 
do not anticipate a need to skip this step for NIRCam 
coronagraphic data. 

During the jump detection step, the jwst pipeline will 
make use of a detection threshold value based on the es- 
timated signal and noise to assess whether a deviation 
between groups is significant enough to be considered a 
jump (see ° for further detail). The default value for 
this threshold is 4, but we repeat an early version of our 
F444W analysis across thresholds of 4 to 16 to search 
for potential improvements. We find that the contrast is 
slightly improved from a threshold of 4 to 5 by ~5x1078 
at 1”, but at larger thresholds it does not vary (devia- 
tions <1x10~® at 1”). Hence, we adopt a detection 
threshold of 5 for all of our NIRCam analyses. 

The Stage 1 products are processed further to Stage 
2 (*calints.fits) files using spaceKLIP, with some addi- 
tional pixel cleaning procedures as follows. Firstly, ev- 
ery pixel with a data quality flag (e.g., indicating hot or 
warm pixels, unreliable data processing) is replaced by 
the median of its orthogonal and diagonal neighbours, 
with the notable exclusion of pixels with a jump flag 
which are typically grouped in clusters. We also inspect 
each pixel for temporal flux variations and identify situ- 
ations where: a) the pixel is bright (MJy/sr > 1) for at 
least one integration, and b) the pixel is relatively faint 
(a >80% decrease in flux compared to the brightest inte- 
gration) for at least one integration. The pixel values for 
the integrations that are not marked as faint are then 
replaced by the median value of the integrations that 
are marked as faint. We note here that an outlier iden- 
tification process based on variations from the standard 
deviation of each pixel in time may be preferable, but 
is difficult to incorporate given the small number of in- 
tegrations in these exposures (see Table 1). Despite the 
above corrections, ~25 static hot pixels remain across all 
images. Although these pixels do not impact our ability 
to recover HIP 65426 b, they introduce residuals in the 
PSF subtraction process (see Section 2.6) that bias our 
measurements of the contrast performance. Therefore, 
we provide the locations of these pixels to spaceKLIP 
manually and correct them in an identical manner as 
the pixels marked with data quality flags. 

As the NIRCam PSF in the F250M filter is under- 
sampled, ringing artifacts are generated by the interpo- 
lation methods in pyKLIP used for image registration 
and spatial shifting of input PSFs. These artifacts bias 


the image registration process, and limit our ability to 
accurately inject synthetic PSFs for contrast curve cali- 
bration and companion astrometry and photometry. To 
overcome this issue, we smooth all of the F250M images 
by a Gaussian filter, as implemented by scipy, with 
o=1.3 pixels. This value of o was chosen as it was the 
minimum possible value that removed the observed arti- 
facts across a test of ten equally spaced values from o=1 
to o=2. We note that this factor of 1.3 is consistent with 
the ratio of the detector pixel scale and the theoretical 
Nyquist sampling at 2.5 ym (assuming a reduced pri- 
mary mirror diameter of 5.2 m due to the NIRCam Lyot 
stop). Smoothing these images may lead to reduced pre- 
cision in our astrometric analysis, however, it should not 
influence the accuracy of our retrieved photometry. 


2.4. MIRI Coronagraphy 


The MIRI observational sequence was executed from 
21:05 July 17th 2022 to 05:19 July 18th 2022 UTC. Ex- 
posures were taken for HIP 65426 b in the F1140C filter 
at one roll angle, then once again at a second roll angle, 
and then for the reference star. This sequence was then 
repeated in reverse for the F1550C filter. This structure 
is slightly different to that of the NIRCam observations 
as each filter is tied to a specific four-quadrant phase 
mask coronagraph (FQPM, Boccaletti et al. 2015), and 
target acquisition must be repeated when switching be- 
tween them. Inserting the reference observations be- 
tween the observations of HIP 65426b minimises the 
time separation between science and reference exposures 
for each filter, and therefore the extent of the wavefront 
evolution between them. After all science and refer- 
ence observations are complete, we perform the dedi- 
cated background observations that are used to subtract 
the dominant “glow stick” stray light feature (Boccaletti 
et al. 2022) as described in Section 2.1. 

We begin data reduction using the Stage 0 (*un- 
cal.fits) files as generated by the jwst pipeline. These 
products are then processed to Stage 1 (*rateints.fits) 
files using spaceKLIP. Similarly to NIRCam (see Sec- 
tion 2.3) we explore the impact of the jump detection 
threshold on our analysis. For these data in particular, 
we found that the default jump detection threshold value 
of 4 is too low and leads to a number of pixels being erro- 
neously flagged as containing a jump. Flagged pixels are 
interpreted differently to unflagged pixels in the ramp 
fitting procedure and the resulting Stage 1 files contain 
a large number of pixels with unrealistic (negative) flux 
values as a result. After repeating an early version of 
our F1140C analysis across thresholds of 4 to 16, we ob- 
serve an improvement in contrast between a threshold 
value of 4 and 5 of ~2x107° at 1”, and ~1x10~° at 2”. 


JWST IMAGING OF HIP 65426B FROM 2—16 MICRON 9 


Beyond this value there is only slight improvement, and 
the obtained contrast varies within ~5x10~° at 1’. For 
our final analyses we select a threshold value of 8, as it 
has the best contrast at 2” and fewer pixels with unre- 
alistic flux values than lower thresholds (as determined 
from visual inspection). 

Following ramp fitting we found that the first integra- 
tion of each exposure contained a significantly increased 
level of noise indicative of a detector reset anomaly. The 
current available calibration dark exposures that are 
used to correct for this anomaly were usually acquired 
in quick succession, and do not have a similar amount 
of dead time (e.g., due to telescope slews/dithers) as 
our science exposures. As a result, during our observa- 
tions the detector electronics were given a much longer 
time to settle, and our integrations exhibit an entirely 
different reset anomaly. This effect is most dominant 
for the first integration of each exposure, but appears 
to persist throughout the entire exposure as well. In 
the future better calibration darks will be provided that 
may resolve this issue, but for the purpose of this study 
we simply exclude the first integration of each exposure 
from all further analysis. 

The Stage 1 products are processed further to Stage 
2 (*calints.fits) files using spaceKLIP, with some addi- 
tional pixel cleaning procedures as follows. First, every 
pixel with a data quality flag (e.g., indicating hot or 
warm pixels, unreliable data processing) is replaced by 
the median of its orthogonal and diagonal neighbours, 
with the notable exclusion of pixels with a jump flag 
which are typically grouped in clusters. Following this 
correction, ~30 static hot pixels remain in our images 
and are corrected following an identical procedure by 
manually providing the pixel locations to spaceKLIP. 
Similarly to NIRCam, these final pixels primarily im- 
pact the measured contrast performance, and not our 
ability to recover HIP 65426 b. 

As shown in Boccaletti et al. (2022), the MIRI corona- 
graphic fields of view are subjected to a stray light “glow 
stick” feature along the horizontal edges of the FQPMs 
that dwarfs the residual stellar flux (see Figure 2). We 
subtract this feature from our processed Stage 2 prod- 
ucts for each filter using a median background image of 
every 4—5 integrations from the dedicated background 
observations to the corresponding 4—5 science or refer- 
ence integrations. The value of 4—5 was selected as it 
provided a slightly improved contrast at 1’/’ compared to 
other tested numbers of integrations per median, rang- 
ing from 1 (~1x1074 improvement) to all available in- 
tegrations (~1x10~° improvement). Grouping the me- 
dian subtraction in this manner better captures the dif- 
fuse reset anomaly noise between integrations mentioned 
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Figure 2. Left: A single integration in the Stage 2 (*cal- 
ints.fits) file for MIRI coronagraphy of HIP 68245 in the 
F1550C filter. Right: As on the left, except following sub- 
traction of a median background frame of an “empty” region 
of the sky. Both images are identically scaled. Before sub- 
traction the residual stellar flux is completely obscured by 
the stray light “glow stick” (Boccaletti et al. 2022), but can 
be easily recovered. 


above. Additionally, it may help capture variations in 
the stray light feature, which varies with a standard de- 
viation of ~0.5% (as identified by variations in the me- 
dian pixel flux for pixels above 50% of the peak pixel 
flux across integrations). 

Following the median background subtraction we find 
that the residual stellar flux is easily recovered (see Fig- 
ure 2). We also attempted to model the stray light 
glow stick as a principal component within both KLIP 
(Soummer et al. 2014) and LOCI (Lafreniére et al. 2007) 
based subtractions, but found that they were suscepti- 
ble to oversubtraction of the residual stellar flux and/or 
could not additionally account for background variations 
between integrations. We anticipate that with careful 
masking and optimisation of these algorithms it may be 
possible to overcome these issues, however, the median 
frame background subtraction is already highly effective 
and improvements to the achieved contrast are likely to 
be minimal. 


2.5. Image Alignment 


The NIRCam and MIRI coronagraphic modes adopt 
independent target acquisition procedures to correctly 
center a star behind each focal plane mask. The in-flight 
positions of the NIRCam mask centers are known to bet- 
ter than ~10 mas but the distortion model is still being 
refined. At present the target acquisition error for the 
MASK335R is as large as ~12—30 mas, or 0.2—0.5 pix- 
els (Girard et al. 2022). This error is dominated by the 
precision of the centering algorithm, which is not well 
adapted to the PSF shape with coronagraphic optics 
(wider in z—axis), and not by the small angle maneu- 
ver (SAM) that places a target at its desired position 
behind the mask (which for NIRCam is repeatable to 
~6 mas). For MIRI the mask center positions have been 
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measured to ~5—10 mas, or ~0.1 pixels (Boccaletti et al. 
2022). However, the SAM has a typical uncertainty of 
~ 10—20 mas, leading to positional offsets between dif- 
ferent rolls or targets. Finally, between integrations the 
pointing stability of JWST (~1 mas) and the accuracy 
of the small-grid dither maneuvers (~2—4 mas) will lead 
to further positional shifts for both the NIRCam and 
MIRI coronagraphs. 

For NIRCam the absolute star position is only ex- 
plicitly measured for the first science image in each 
filter. This position is measured using a cross cor- 
relation of a model coronagraphic PSF as obtained 
from webbpsf_ext to the science PSF using the 
scikit-image package (van der Walt et al. 2014). To 
best match the science PSF, the model PSF is generated 
using the telescope optical path difference (OPD) map 
as obtained on July 29th. To identify the accuracy of 
this process we repeat the procedure, except comparing 
the model PSF to itself, and to a second model PSF that 
was generated using a pre-launch measurement of the 
telescope OPD (which exhibits comparable differences 
to the model PSF as our data). In each case, we man- 
ually shift the comparison PSFs across a range of 0.01 
to 0.5 pixels and attempt to recover these offsets using 
the cross correlation process. For the self comparison, 
we recover the injected shift to at least ~0.03 pixels, or 
1.9 mas, whereas for the second model PSF comparison, 
we recover the injected shift to ~0.1 pixels, or 6.3 mas. 
Given the comparable difference between the model PSF 
for both the second model PSF and our data, we adopt 
the latter uncertainty as a systematic uncertainty in our 
astrometric measurements (see 4.3). The shifts of the 
other science and reference images are obtained relative 
to the first science image through a similar cross cor- 
relation procedure. However, we instead use a box of 
11x11 pixels around the coronagraph center position, 
where the flux is dominated by the central core of the 
coronagraphic PSF. All measured relative shifts match 
expectations for the JWST pointing precision of ~1 mas 
(lo, radial) (Rigby et al. 2022). 

To estimate the absolute star position for MIRI, 
we use the current measurements for the centers of 
the coronagraphic masks and assume the the star 
is perfectly centered behind the coronagraphic mask. 
In zero-indexed subarray x-y coordinates these val- 
ues are (119.749, 112.236) and (119.746, 113.289) for 
the FQPM1140 and FQPM1550, respectively (Jonathan 
Aguilar, private communication). As a result, the stel- 
lar position is not known better than a minimum of 
~10 mas (Boccaletti et al. 2022). Similarly, the rel- 
ative alignment between images may be discrepant by 
~1—4 mas (Rigby et al. 2022), or <0.05 pixels. At- 


tempts to estimate the absolute star position in a simi- 
lar manner to the NIRCam data were unsuccessful, and 
this process led to significantly larger estimated shifts 
than the known pointing stability of JWST (~50 mas, 
vs ~1 mas). This is most likely due to the variable 
spatial structure of the residual PSF, which signifi- 
cantly changes with small pointing shifts. An effort 
was also made to fit each MIRI image with model PSF's 
across all feasible pointings using webbpsf_ext® within 
an emcee (Foreman-Mackey et al. 2013) Markov chain 
Monte Carlo (MCMC) framework. However, these fits 
were unable to converge and upon visual inspection the 
spatial structure of the empirical and model PSFs were 
different, despite modelling the PSFs based on measure- 
ments of the telescope OPD within ~1 day of our obser- 
vations. We do not believe that using the coronagraphic 
mask centers as a proxy for the absolute star position 
has significantly affected our results, but it is certainly 
an area of improvement for future studies using MIRI 
coronagraphy. 

For NIRCam all images are aligned to a common cen- 
ter based on the measured shifts, however, for MIRI 
we opt to not perform any realignment. This decision 
is made under the assumption that because a pointing 
shift does not primarily cause a translation of the resid- 
ual PSF in the MIRI images (in contrast to NIRCam), 
the unshifted reference images are more descriptive of 
variations in the science images. When comparing a 
RDI subtraction using unshifted reference images, and 
a separate RDI subtraction with a realignment based 
on the ideal small-grid dither positions, the measured 
contrasts are in agreement and this choice does not sig- 
nificantly impact our results. 


2.6. PSF Subtraction 


We perform a subtraction of the residual stellar PSF 
in each filter following three different principal compo- 
nent analysis (PCA) based methods as implemented in 
spaceKLIP. First, we take the two independent rolls of 
HIP 65426 b and perform an ADI subtraction. Second, 
we perform an RDI subtraction by using the correspond- 
ing observations of the reference star, HIP 68245, as a 
PSF library. Finally, we perform an ADI+RDI subtrac- 
tion, which is identical to the RDI subtraction except 
that images at the opposite roll angle are also included 
in the PSF library. In each case the subtraction is per- 
formed on each integration from both science rolls indi- 
vidually, before being rotated to a common orientation 
as marked in Figure 3 and summed together. Although 


8 https://github.com/JarronL/webbpsf_ext 
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Figure 3. Unsubtracted and KLIP subtracted image stamps for the NIRCam F356W (top row) and MIRI F1140C (bottom 
row) filters. The leftmost column displays the median unsubtracted image for a single science roll, and all other columns display 
the KLIP subtracted images for ADI, RDI, and ADI+RDI subtraction methods using the maximum number of KLIP PCA 
modes. All images are oriented as shown by the directional arrow in unsubtracted image column, the position of the star (white 
star) is also marked. Additionally, the intensity of all images for a given filter are identically scaled. The exoplanet, HIP 65426 b, 
can be easily identified at a position angle of ~150° in the subtracted images. We note that the distinct “hamburger” shaped 
central core and six-lobed structure of the companion PSF in the NIRCam images is an expected feature that is related to the 
Lyot stop design, and not indicative of discrete astrophysical sources. 


the number of annuli and subsections the PSF subtrac- 
tion is performed across can be adjusted, we find that 
this does not improve the observed contrast. Hence, we 
perform all subtractions using a single annulus and a 
single subsection (i.e., the entire image). The number of 
KLIP PCA modes can also be adjusted to tune the ag- 
gressiveness of the PSF subtraction. Hence, we perform 
the PSF subtraction across the full range of possible 
PCA modes to investigate the impact on our measured 
contrast and companion fitting. The maximum number 
of PCA modes is dependent on the exposure settings 
for each filter and corresponds to: the number of inte- 
grations in a single roll for ADI, the total number of 
integrations across all 9 dithers for RDI, and the sum of 
the two for ADI+RDI (see Table 1 for precise values). 

Pre- and post-subtraction images for the NIRCam 
F356W and MIRI F1140C filters are shown in Figure 3, 
and images for all filters are shown in Appendix A. We 
note that the distinct “hamburger” shaped central core 
and six-lobed structure of the companion PSF in the 
NIRCam images is an expected feature that is related 
to the Lyot stop design, and not indicative of discrete 
astrophysical sources. 


3. ANALYSIS 
3.1. Contrast Calibration 


All proceeding contrast measurements are determined 
relative to a synthetic spectrum of HIP 65426 in each of 
the JWST filters, as estimated from fitting stellar and 
disk models to existing photometry following Yelverton 
et al. (2019), see Figure 4. We use data of HIP 65426 
from Hipparcos/Tycho-2 (Hg et al. 2000), Gaia DR2 
(Gaia Collaboration 2018), 2MASS (Cutri et al. 2003), 
ALLWISE (Wright et al. 2010), AKARI IRC (Ishihara 
et al. 2010), and Spitzer MIPS (Chen et al. 2012). 
The fitting procedure compares synthetic photometry 
of models to the data to compute a x? value, and pos- 
terior distributions are found using MultiNest (Feroz 
et al. 2009; Buchner et al. 2014). We derive our own 
zero points using the CALSPEC Vega spectrum (Bohlin 
et al. 2014). We use PHOENIX models (Allard et al. 
2012) for the stellar photosphere, and a Planck func- 
tion for the disk model. There is a small excess at 
24 ym that was previously reported at 3.50 by Chen 
et al. (2012), though not considered significant in that 
paper. The best fit model has an effective tempera- 
ture of 8600+200 K and luminosity 16+1 Lo. The dust 
temperature and luminosity are very poorly constrained 
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Figure 4. The best fit stellar model (purple) to existing: 
Hipparcos/Tycho-2, Gaia DR2, 2MASS, ALLWISE, AKARI 
IRC, and Spitzer MIPS data for HIP 65426 (black circles). 
The stellar (blue) and disk (red) components of the model 
are also shown, along with the equivalent model fluxes in 
each photometric band (solid circles). Error bars are smaller 
than the circle diameter, except for the 70 um MIPS point, 
which is instead marked as an upper limit (black triangle). 


(Taust = 3007309 K, and Laust/L. ~ 2 x 10-5), though 
this uncertainty does not significantly impact the flux 
estimation in JWST bandpasses because the excess is 
small. The flux excess at 24 wm is not high, but if real 
the dust would reside relatively close to the star, prob- 
ably less than 1 au (hence we use the total model flux 
to compute the contrast, as any dust component that 
contributes IR flux would remain unresolved). We use 
the posterior distribution of model parameters and syn- 
thetic photometry to generate a distribution of fluxes in 
the JWST bands, and adopt the maximum likelihood 
solution for the stellar flux in each bandpass. 


3.2. Contrast Curves 


Following PSF subtraction, we obtain metrics of 
the sensitivity as a function of angular separation 
(i.e., contrast curves) for all observational filters us- 
ing spaceKLIP. To avoid biasing the contrast measure- 
ment we mask regions of the subtracted images near 
HIP 65426 b, background sources, and the FQPM edges. 
These “5c” contrast curves report the flux level corre- 
sponding to a 50-equivalent false alarm probability of 
2.9 x 10~7 after correcting for small sample statistics at 
small separations (Mawet et al. 2014). We call them 5c 
contrast curves for brevity. 

To obtain a more accurate measurement of the con- 
trast performance, the throughput of the coronagraph 
and the intrinsic throughput of the KLIP subtraction 
must be accounted for by injecting and then recovering 
the flux of artificial sources (Adams & Wang 2020). All 


artificial sources are generated using webbpsf? (Perrin 
et al. 2012, 2014) at an initial intensity equivalent to 
a signal-to-noise ratio of 25. Immediately prior to in- 
jection and based on a desired injection location, each 
source is modulated by the coronagraphic throughput 
using a synthetic throughput map. 

Synthetic coronagraphic throughput maps are pro- 
vided in the calibration reference files!?, however, both 
the provided NIRCam and MIRI FQPM maps are in- 
accurate. For the NIRCam MASK335R, the position 
of the occulting masks within the throughput map does 
not correspond with its actual location. Therefore, we 
modify the throughput map by extracting the pixels im- 
pacted by the occulting mask and repositioning them at 
the true mask center location of (149.9, 174.4) in zero- 
indexed subarray x-y coordinates (Jarron Leisenring, 
private communication). In the case of MIRI, all of the 
FQPM maps are rudimentary and do not accurately cap- 
ture the spatial throughput variations. Therefore, we in- 
stead use custom simulated maps of the FQPM through- 
put produced using webbpsf_ext. In brief, webbpsf_ext 
generates position-dependent PSF'’s based upon interpo- 
lation over a densely-sampled grid of individual PSFs 
along the FQPM axes. A uniform (flat) illumination 
of 10-7 Jy was convolved with the position-dependent 
PSFs to produce an illumination pattern that was then 
normalised and scaled to the subarray dimensions for 
each MIRI FQPM. 

Once scaled by the coronagraphic throughput, PSF's 
are injected into multiple copies of the unsubtracted 
science images across a range of separations extend- 
ing to 4”, and for a range of position angles from 0 
to 360°. Sources are not injected within 2A/D of each 
other or a masked region. These images then undergo 
KLIP subtraction in an identical manner to the sci- 
ence images, and the relative flux of an initial source 
PSF and the KLIP processed PSF as estimated within 
pyKLIP describes the overall coronagraphic mask plus 
KLIP throughput at each location. Finally, the basic 
50 contrast is divided by an interpolation of the me- 
dian throughput across all position angles to obtain the 
calibrated contrast. 

For the ADI subtraction, the measured contrast does 
not vary significantly when using more than 1 PCA 
mode for the NIRCam filters and more than 2 PCA 
modes for the MIRI filters. For RDI, the measured con- 
trast does not vary beyond ~6 modes for both NIRCam 
and MIRI. Finally, for ADI+-RDI we find that there is a 


9 https: //webbpsf.readthedocs.io/en/latest / 
10 https: //jwst-crds.stsci.edu 
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transition to improved contrast at Pnax — Pap1 modes, 
where Pyax is the maximum possible number of PCA 
modes possible, and Pap, is the maximum number of 
modes for the ADI subtraction. This is likely a result of 
the much larger number of reference images weighting 
the calculation of the PCA modes to be mostly RDI-like, 
until a sufficient limit is reached where the influence of 
the opposing roll images appears in the principal com- 
ponents. Beyond this transition value, the measured 
contrast does not vary significantly. 

Example calibrated contrast curves for the NIRCam 
F356W and MIRI F1140C filters using the maximum 
number of PCA modes are shown in Figure 5, and con- 
trast curves for all filters are shown in Appendix B. 


3.3. HIP 65426b PSF Fitting 


To analyse the properties of HIP 65426b in greater 
detail, we make use of the forward model PSF fit- 
ting routine provided by pyKLIP and implemented in 
spaceKLIP. Briefly, this routine takes a model of the 
companion PSF and uses a forward model of the KLIP 
subtraction to apply PSF distortions that arise natu- 
rally from the KLIP process. This resultant PSF can 
then be scaled/shifted to best match the observed com- 
panion PSF and obtain a measurement of its location 
and intensity (Pueyo 2016; Wang et al. 2016). 

For our analysis we adopt an independent model PSF 
for each filter using webbpsf_ext functionality as imple- 
mented in spaceKLIP. Specifically, we use webbpsf_ext 
to generate an offset PSF at the predicted location of 
HIP 65426 b as adopted from whereistheplanet (Wang 
et al. 2021) and at the appropriate position for each sci- 
ence roll image. Each NIRCam and MIRI PSF is gen- 
erated using a measured OPD map as determined from 
wavefront sensing and control observations on July 29th 
2022 and July 17th 2022, respectively. Particularly for 
MIRI, the PSF of an off-axis source is sensitive to its 
location relative to the FQPM edges, and explicitly gen- 
erating a model PSF close to this location is important 
for obtaining a close match to the true companion PSF. 
As the spatial intensity of the model PSF depends on an 
assumed spectral energy distribution (SED), we use an 
existing best fit BT-SETTL model for HIP 65426b from 
Cheetham et al. (2019). We note that as our initial 
model PSF is normalised to unit intensity, the primary 
purpose of selecting a model SED is to capture the rela- 
tive variation in flux as a function of wavelength across 
each bandpass, and is not meant to accurately predict 
the absolute flux of HIP 65426 b in each bandpass. 

The input model PSF is converted to physical units by 
scaling it to the flux of HIP 65426 as estimated by the 
stellar model described in Section 3.1. Therefore, all 


derived photometry is anchored relative to our assump- 
tion of the stellar flux. Furthermore, any comparisons 
between the intensity of this PSF and the observed PSF 
is subject to an implicit assumption that the absolute 
flux calibration of JWST is perfect. In reality, the abso- 
lute flux calibration accuracy requirement for NIRCam 
and MIRI coronagraphy is 5% and 15%, respectively!'. 
Both these fundamental uncertainties on the absolute 
flux calibration and the uncertainty on the stellar model 
flux as derived from the distributions described in Sec- 
tion 3.1 are propagated as an increased uncertainty on 
all of our flux measurements. 

For each filter we fit the location and intensity of the 
model PSF to the true PSF of HIP 65426b from the 
ADI+RDI subtraction using 6 PCA modes. The fitting 
procedure is executed in an MCMC framework as imple- 
mented by emcee (Foreman-Mackey et al. 2013), with 
50 walkers for 300 steps, of which the first 100 steps are 
discarded as burn-in. spaceKLIP allows for the spatial 
scale of noise in an image to be fit with a variety of 
different kernels as implemented in a Gaussian process 
framework. For NIRCam the noise in our images ap- 
pears to be correlated at the separation of HIP 65426 b, 
so, following Wang et al. (2016), we initially adopt a 
Matérn 3/2 kernel to better capture this spatially cor- 
related noise structure. However, the generated model 
PSF is slightly mismatched to the observed data and a 
positive flux residual was present at the PSF core. Fu- 
ture improvements to model PSF generation may allevi- 
ate this issue, but in this work we instead assume uncor- 
related noise (using the “diagonal kernel” option), which 
is able to better capture the flux of the companion at 
the expense of underestimating the obtained error bars. 
Therefore, more realistic error bars for the NIRCam pho- 
tometric measurements are determined through a pro- 
cess of companion injection and recovery. For each fil- 
ter, the best fit model PSF is used to subtract away 
the companion flux, and is then injected at 20 different 
position angles spanning 0—105° and 195—360° across 
an equivalent number of duplicate science images. The 
HIP 65426b PSF fitting procedure is then repeated on 
these synthetic PSF's and the standard deviation in the 
measured flux across all 20 position angles is adopted as 
the estimated error. In the case of MIRI the observed 
noise is visually consistent with uncorrelated noise, so 
we adopt a diagonal kernel here as well. 

Results from the PSF fitting procedure are discussed 
further in Section 4, and images of the data, model, and 
residuals to each PSF fit are displayed in Appendix C. 


11 https: //jwst-docs.stsci.edu/jwst-data-calibration-considerations 
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Figure 5. Contrast curves for observations in the F356W and F1140C filters using both an ADI (dotted lines), RDI (dashed 
lines) and ADI+RDI (solid lines) subtraction with the maximum possible number of PCA modes. Both the measured contrast 
of the true on-sky observations (black lines) and predicted contrasts as generated from PanCAKE (Carter et al. 2021b) (light 
blue lines) are displayed. The effective inner working angles (IWA), corresponding to the separation at 50% transmission, are 
also marked (red dashed lines). Contrast curves for all other filters are displayed in Figure 14. 


4. DISCUSSION 


The known companion HIP 65426 b is clearly detected 
in all seven of the observational filters using RDI, and 
all filters except the MIRI F1550C using ADI. These ob- 
servations represent the first images of an exoplanet to 
be obtained with JWST, and the first ever direct detec- 
tion of an exoplanet beyond 5 wm. Furthermore, these 
observations provide a first look at the achieved contrast 
for JWST’s high contrast imaging modes following the 
completion of observatory commissioning. 


4.1. Achieved Contrast 


These observations demonstrate that, in addition to 
the existing work of Girard et al. (2022), Kammerer 
et al. (2022), and Boccaletti et al. (2022), both NIR- 
Cam and MIRI coronagraphic imaging are exceeding 
their predicted contrast performance. Examples of this 
improvement in comparison to pre-launch contrast es- 
timates as obtained by PanCAKE (Carter et al. 2021b) 
are demonstrated in Figure 5, and contrast curves for 
all seven filters are displayed in Appendxi B. 

At the MASK335R nominal IWA of 0.64”, the F356W 
data achieve a contrast of ~2x107°, sloping down to 
~4x107® at 1” and then ~9x10-7 beyond 3”. In com- 
parison, at the FQPM1140C nominal IWA of 0.36”, the 
MIRI F1140C data achieve a contrast ~1x1072, slop- 
ing down to ~2x10~4 at 1” and then ~5x10~° beyond 
3”. For brevity we do not describe the achieved con- 


trast in the other filters, and instead refer the reader to 
Appendix B. 

In the background limited regime beyond ~2” our 
measured contrast approximately matches the predicted 
sensitivity for NIRCam ADI and ADI+RDI subtrac- 
tions, and is ~2 times more sensitive than the prediction 
for the RDI subtraction. For MIRI in this regime, all 
three subtraction methods outperform their predicted 
sensitivity by a factor of 1.5 — 2. In the contrast lim- 
ited regime below ~2” we observe similar improvements 
upon the predicted contrast, with the ADI and RDI sub- 
tractions demonstrating a factor of up to 10 times deeper 
contrast, and the ADI+RDI subtraction improving by 
a factor of ~2. In some instances at the shortest sepa- 
rations below ~0.4” the contrast does underperform by 
up to a factor of ~2 compared to predictions for RDI 
subtractions. The primary driver for these improve- 
ments is likely the improvement in the overall optical 
and stability performance of JWST compared to expec- 
tations. The total throughput is ~10—20% larger than 
predictions, driving analogous improvements in signal- 
to-noise; the overall telescope wavefront error is ~ 50% 
smaller than requirements (75/110 nm vs 150/200 nm 
for NIRCam/MIRI), improving the raw contrast; and 
the pointing stability of 1 mas is ~6—7 times smaller 
than predictions, meaning smaller drifts in alignment 
behind the coronagraphic mask throughout an observa- 
tion (Rigby et al. 2022). 
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The different contrasts as achieved by the ADI, RDI, 
and ADI+RDI subtractions allow for more concrete rec- 
ommendations in observing structure for future pro- 
grams. Most significantly, the improvement between 
RDI and ADI+RDI subtractions is negligible. Future 
observers will may be able to achieve their science goals 
with less telescope time by focusing on purely ADI or 
RDI subtraction strategies, however, a broader range 
of observations will be required to rule out ADI+RDI 
strategies entirely. Within 1—2” the ADI subtraction for 
MIRI F1550C, and to some extent F1140C, struggles to 
fully subtract the residual stellar PSF and an RDI based 
subtraction strategy should be preferred. Although ADI 
subtractions provide demonstrably worse contrast across 
all filters, they require a significantly smaller amount 
of observing time due to the cost-intensive nature of 
dithered reference observations. However, similar con- 
trasts to ADI can be achieved with an undithered ref- 
erence observation as well (Girard et al. 2022). To re- 
duce the overall telescope time required for high contrast 
imaging programs, future observers should carefully con- 
sider the contrast performance necessary to meet their 
science goals, and adjust their PSF subtraction strategy 
accordingly. 

Beyond this work, there is significant potential for the 
contrast performance for both NIRCam and MIRI to 
improve. Kammerer et al. (2022) have already demon- 
strated that the NIRCam bar masks can provide deeper 
contrasts at shorter separations, if a full 360° field-of- 
view is not required. Similarly, in later cycles it may 
be possible to position a star at the “NARROW?” bar 
mask position in an attempt to reduce the effective in- 
ner working angle (IWA) even further. The efficacy of 
RDI subtractions may improve with the use of a larger 
PSF library populated by on-sky observations across 
multiple programs. It may also be possible to perform 
an effective PSF subtraction using a much larger PSF 
library composed either entirely or partially of model 
PSFs. A particular opportunity for improvement will 
also come from JWST program GO-02627 (Ygouf et al. 
2021), which aims to estimate the on-sky instrumental 
aberrations that drive variations in the observed PSF 
structure with a model based phase retrieval algorithm. 


4.2. Mass Sensitivity 


Using the obtained contrast curves we also determine 
the detectable mass limits for our observations following 
the approach of Carter et al. (2021la). Briefly, we con- 
vert our contrast curves to magnitude sensitivity curves 
using the stellar magnitudes as described in Section 3.1, 
and then convert these to a mass sensitivity using an 
interpolation of the evolutionary models of Linder et al. 


(2019) (BEX) and Phillips et al. (2020) (ATMO) as- 
suming an age of 14+4 Myr (Chauvin et al. 2017) and 
distance of 107.49+0.40 pce (Gaia Collaboration 2020). 
As in Carter et al. (2021a) we select the chemical equi- 
librium, non-cloudy models to maintain model consis- 
tency across mass ranges. Clouds and disequilibrium 
chemistry likely play a significant role in sculpting the 
emission of sub-stellar atmospheres, however, an inves- 
tigation into these effects is beyond the scope of this 
work. 

Following the calculation of mass sensitivity curves we 
use the Exoplanet Detection Map Calculator (Exo-DMC, 
Bonavita et al. 2012, 2013; Bonavita 2020)!? to estimate 
detection probability maps. In this case, we produce 
a population of synthetic companions with masses and 
semi-major axes from 0.1M)j,, to 100Mjyp and 0.lau 
to 10,000 au respectively. All orbital parameters are 
uniformly distributed except for the eccentricity, which 
is generated using a Gaussian distribution with w = 0 
and o = 0.3 (excluding negative eccentricities; Hogg 
et al. 2010). Implicit in this is an assumption that these 
synthetic companions do not necessarily have a similar 
inclination to HIP 65426 b. The resulting map takes into 
account the effects of projection when estimating the 
detection probability, and the probability for a potential 
companion to truly lie in the instrumental field of view. 

For these particular observations of HIP 65426, we 
identify the F444W and F1140C filters as the most sen- 
sitive to the lowest mass companions for NIRCam and 
MIRI, respectively, and display their detection probabil- 
ity maps in Figure 6. Of the two, the F444W filter is the 
most sensitive, reaching a sub-Jupiter mass sensitivity 
from ~150—2000 au at a 50% probability, and a mini- 
mum sensitivity of ~0.4 Mj,, from ~300—1500 au at a 
50% probability. In contrast, the F1140C is unable to 
reach sub-Jupiter mass sensitivity, and has a minimum 
sensitivity of 1-2 Mjyp from ~150—2000 au. As we 
detect no sources within the observed field-of-view that 
have colours consistent with a planetary mass compan- 
ion, we set equivalent limits on the presence of additional 
companions in the HIP 65426 system. 

As discussed in Carter et al. (2021a), at a given dis- 
tance A stars are generally poor targets for detecting 
the lowest mass planets in terms of detection sensitiv- 
ity, whereas M stars are among the best. Given the im- 
proved performance of JWST, it is likely that for targets 
within (or outside of) the M star sample from Carter 
et al. (2021a) it will be possible to detect Uranus and 
Neptune mass objects beyond ~100—200 au, and Sat- 


12 https://github.com/mbonav/Exo_DMC 
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Figure 6. Detection probability maps as generated by 
Exo-DMC for the most sensitive NIRCam (F444W, top) and 
MIRI (F1140C, bottom) filters. Solid black contours signify 
the 10%, 50%, 85%, and 95% detection thresholds. 


urn mass objects beyond ~10 au. Initial searches across 
a small sample of stars for sub-Jupiter mass objects will 
be performed as part of guaranteed time observations 
(Schlieder et al. 2017). Furthermore, for the nearest tar- 
gets we may be able to push these limits even further, 
with planned observations of a Cen A aiming to be sen- 
sitive to 5 Rg companions from 0.5—2.5 au (Beichman 
et al. 2020, 2021). 


4.3. Astrometry 


The measured astrometry in each of the observed fil- 
ters is obtained from the ADI+RDI reduction and shown 


in Table 2. Each of the measured uncertainties for the 
NIRCam astrometry are propagated with an additional 
uncertainty of 6.3 mas (0.1 pixels), and the uncertain- 
ties of the MIRI astrometry are propagated with an ad- 
ditional uncertainty of 10 mas (0.1 pixels) to account 
for the assumed precision of the absolute star centering 
as described in 2.5. All NIRCam and MIRI astrometry 
are consistent within lo, and in combination the NIR- 
Cam and MIRI astrometry provide a measurement of 
the separation, p=820+6 mas, and the position angle, 
6=149.5+0.4°. To compute these average values, we do 
not treat the NIRCam filters as independent and instead 
average both the quantities and their uncertainties. The 
absolute position of the star on the detector does not 
change significantly (<1 pixel) between NIRCam filters, 
and it is feasible that the measured position has a simi- 
lar systematic offset in each filter (see Section 2.5). As 
the dominant noise source for the NIRCam alignment is 
this systematic offset, it is not appropriate to propagate 
the uncertainties of the NIRCam astrometry in a typical 
fashion. 

We combine these new measurements with the ex- 
isting astrometry from Chauvin et al. (2017) and 
Cheetham et al. (2019) to determine updated orbital 
parameters using the orbitize! package (Blunt et al. 
2020). As in Cheetham et al. (2019), we exclude 
the NaCo epochs due to the inconsistency with the 
SPHERE epochs. Additionally, we exclude the MIRI 
epoch as the observations were taken just two weeks af- 
ter the NIRCam observations and have larger uncertain- 
ties. orbitize! is initialised assuming one companion 
to the primary (HIP 65426b), a total system mass of 
1.97£0.046 Me (Chauvin et al. 2004), and a parallax 
of 9.303140.0346 mas (Gaia Collaboration 2020). Orbit 
generation is performed using the Orbits for the Impa- 
tient (OFTT) algorithm (Blunt et al. 2020) until 100,000 
possible orbits are identified. A random sample of 100 
of the possible orbits along with posterior distributions 
for the entire sample are shown in Figure 7. We are 
able to constrain the semi-major axis, a = Cale i au, 
and the inclination, i = 100+, degrees. Additionally, 
as the motion of HIP 65426b is primarily radial, solu- 
tions that place the line of nodes close to its position 
angle are preferred, and the position angle of nodes is 
also constrained to two possible solutions. 

The addition of the NIRCam astrometry does 
not significantly improve the orbital constraints for 
HIP 65426 b, and all retrieved properties are consistent 
with those from Cheetham et al. (2019) and Bowler 
et al. (2020). Although our eccentricity distribution 
more strongly favours higher eccentricities compared to 
Cheetham et al. (2019), it remains essentially uncon- 
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Figure 7. Orbital fitting of both JWST NIRCam (this work) and SPHERE (Chauvin et al. 2017; Cheetham et al. 2019) 
astrometric measurements of HIP 65426b using orbitize! (Blunt et al. 2017). Top Left: A random sample of 100 orbit 
models from the retrieved posterior (purple curves). The position of the planet (white circle) and the star (white star) are 
marked, and the epoch at a given position in an orbit is indicated by the colour bar. Top Right: Separation and position 
angle versus epoch for both the SPHERE (squares) and JWST NIRCam astrometry (circles). The same random sample of 100 
orbits is also displayed (grey lines). Bottom: Posterior distributions for the semi-major axis (a), eccentricity (e), inclination 
(i), argument of periastron (w), and position angle of nodes (Q). The 50% quantile from these distributions (green dot-dashed 
line) are also indicated. In this particular case, these additional JWST observations do not significantly increase the constraints 
on the measured orbital parameters. 


strained and should not be interpreted as evidence for Collaboration et al. 2019; Lacour et al. 2021; Hinkley 
a highly eccentric orbit for HIP 65426b. However, if et al. 2022b) and has an even greater potential to im- 
this high eccentricity is real, it would give credence to prove upon our reported constraints. 


the scenario proposed in Marleau et al. (2019), where 
HIP 65426 b initially formed through core accretion be- 


: F : ; 4.4. Photomet 
fore being scattered to a wider separation by an addi- aia 


tional companion. The ability of JWST to provide high As with the astrometry, the measured photometry in 
precision astrometry may improve with improvements all of the observed filters is obtained from the ADI+RDI 
to the measurement of the absolute star position in the reduction and shown in Table 2, the subtracted images 
NIRCam images, which in this case is limited by the pre- of all of these filters is shown in Figure 8, and the pho- 
cision with which we can locate the star center behind tometric data points themselves are shown in Figure 9. 
the coronagraphic mask (see Section 2.5). The measured contrast of the planet relative to the star 

Separate to JWST, HIP 65426 b has been observed as ranges from 10.403 mag in the F250M filter to 8.445 mag 
part of the ExoGRAVITY program (Lacour et al. 2020, in the F1550C filter. Images for the ADI and RDI sub- 
Program ID: 1104.C-0651), which has routinely demon- tractions can be found in Appendix A. Additional liter- 
strated sub—mas astrometric precision (e.g., Gravity ature photometric measurements as shown in Figure 9 


are provided in Appendix D. 
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Figure 8. Images of the exoplanetary companion, HIP 65426 b, in all seven NIRCam and MIRI filters used in our observations. 
Each image is produced following an ADI+RDI KLIP subtraction of the residual stellar PSF. The measured position of the star 
is marked (white stars), and the orientation and pixel scales of all images are marked in the top left panel. 


Filter —_p (mas) 6 (deg) m. (mag) 


A (mag) 


F250M 820£7 150.140.5 6.78340.054 10.403+0.046 
F300M 817£7 = =149.7+0.5 6.7660.046 10.014+0.043 
F356W 80847 149.7+40.5 6.767+0.048 9.280+0.021 
F410M 8137 149.7+0.5 6.765+0.051 8.978+0.020 


F444W 9 815+7 = 149.8+0.5 6.764+0.054 8.939+0.019 
F1140C =827+11 149+1 6.722+0.038  8.653+0.018 
F1550C 828+15 149-1 6.7660.072  8.445+0.035 


Acorr (mag) mp (mag) Flux (Wm7? m7 ) 
10.403+0.071  17.186+0.088 (3.3540.27) x 107 is 
10.014+0.069 16.780+0.083 (2.44+0.19) x 107 ii 
9.28040.058  16.04740.076  (2.55+0.18) x 1072” 
8.978£0.058  15.743+0.077  (1.99+0.14)x 10717 
8.9390.058  15.703+0.078  (1.59+0.11)x 107!” 
8.653£0.164  15.375+40.168  (5.17-40.80)x 10-9 
8.445+0.167 15.211+0.182 (1.72+0.29) x 107 7 


Table 2. JWST astrometry and photometry of HIP 65426b, m. corresponds to the stellar magnitude in each filter, and 
Acorr corresponds to the relative magnitude following the propagation in uncertainties of a 5% or 15% absolute flux calibration 
accuracy for NIRCam and MIRI, respectively. The position angle (0) is provided from North through East, and all apparent 


magnitudes are relative to Vega. 


A notable tension can be seen between the JWST 
NIRCam data shown here and the VET NaCo L’, 
NB4.05, and M’ data (Cheetham et al. 2019; Stolker 
et al. 2020b, Appendix D). The model fit shown in Fig- 
ure 9 (also see Section 4.7) does not include the NaCo 
data, and matches the JWST NIRCam data well. The 
VLT NaCo L’ and M’ photometric points are the aver- 
age of the values reported in Cheetham et al. 2019 and 
Stolker et al. 2020b, and are discrepant with the predic- 
tions of this model by 2.10 and 2.70, respectively. In 
light of this, and given the recent advent of high con- 
trast imaging with JWST, we explore the potential for 
either data set to be biased in some way. 


In the case of JWST, we focus our efforts on inves- 
tigating the flux calibration of NIRCam coronagraphy. 
The required absolute flux prediction accuracy for this 
mode is 5%!%, but has yet to be formally verified from 
calibration observations. First, we perform simple aper- 
ture photometry on calibration observations of the star 
2MASS J18022716+6043356 (A2V, 2MASS K,=11.83) 
from CAL-01536 (PI: K. Gordon) with NIRCam F356W 
imaging. The retrieved flux is ~7% higher than a 
synthetic F356W observation of an A2V PHOENIX 


13 https: //jwst-docs.stsci.edu/jwst-data-calibration-considerations 
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Figure 9. All existing spectroscopic and photometric observations of HIP 65426 b as obtained from SPHERE/IFS (triangles), 
SPHERE/IRDIS (squares), NaCo (diamonds), and JWST (circles). Top: Data are plotted alongside the 1, 2, and 30 confidence 
intervals obtained from fitting to a collection of BT-SETTL atmospheric forward models (blue shaded regions), and the model 
values in the photometric bandpasses (small blue circles). At 30, the best fit models occupy parameter ranges of Teg = 167372! K, 
log(g) = 410 dex, and R = 0.9073 04 Rjup. The NaCo data have not been included in the model fitting process. Also 
plotted are the normalised filter throughput profiles for all photometric observations, with the NaCo throughputs scaled by a 
factor of 2 to improve clarity. Bottom: Residuals of each data point relative to the best fit model in addition to 1, 2, and 30 


regions (grey shading). 


model (Husser et al. 2013) which was normalised to the 
2MASS kK, magnitude listed above. Secondly, follow- 
ing a similar approach we compare the flux of sources 
that were observed in both NIRCam Imaging and NIR- 
Cam Coronagraphy from COM-01070 (PI: J. Girard) in 
the F335M filter. We selected three sources from the 
coronagraphic observations that were: underneath the 
mask substrate, not near an occulting mask, not over- 
lapping with a nearby source, and present in the imag- 
ing observation field-of-view. The fluxes as estimated 
from the coronagraphic observation were 1%, 4%, and 
7% higher than those estimated from the imaging ob- 
servation. Given the simplicity of this analysis, we do 
not attempt to precisely quantify the accuracy of the 
current flux calibration, however, these measurements 
provide strong evidence that it is not driving the ob- 


served tension with the NaCo data, which is a factor of 
~2 larger than the analogous JWST data. 

With reference to the NaCo data, using indepen- 
dent analyses Cheetham et al. (2019) and Stolker et al. 
(2020b) extract photometry consistent with each other 
to lo. Both the LZ’ and M’ observations were taken 
in average weather conditions, with light but varying 
cloud cover!*. The non-contemporaneous nature of flux 
calibrations for ground-based observations can lead to 
small biases if the weather conditions change between 
flux and science frames. Further, the M’ observation is 
subject to a high thermal background, and the residual 
intensities near the location of the companion appear 
to be approximately half of the companion peak inten- 


14 https: //www.eso.org/asm/ui/publicLog 
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sity (Stolker et al. 2020b), though the impact of the 
high background is seen in the error bar of the M’ pho- 
tometry. The NB4.05 observation was taken in particu- 
larly poor and variable conditions, with a mean seeing of 
1.52+3.01” (Stolker et al. 2020b), hence also susceptible 
to flux variations. While the points raised above sug- 
gest it is conceivable for the NaCo photometry to have 
been overestimated, pinpointing the exact origin of the 
discrepancy would involve an in-depth re-analysis of the 
NaCo data, which is left for future work. 

The observed tension may also be influenced by gen- 
uine variability of HIP 65426b, its host star, or both. 
However, we do not investigate this such a scenario in 
this work. 


4.5. Bolometric Luminosity 


With the addition of JWST NIRCam and MIRI pho- 
tometric observations, the SED of HIP 65426b is mea- 
sured over a significant portion of the range from 1 wm 
to 15 um. These measurements enable a tight constraint 
on the bolometric luminosity of the planet. 

To calculate the luminosity, a full SED was created 
by distributing the flux-density from photometric mea- 
surements over the effective bandwidth for each filter 
and using a model atmosphere to extrapolate beyond 
and interpolate between measured bands. Luminosity is 
then determined by integrating this semi-empirical SED 
over wavelength. 

Since all of the flux measured in the NIRCam/F410M 
photometry is accounted for in the F444W measure- 
ment, we used only the wider band for our analysis. We 
also added measurements from the literature at shorter 
wavelengths, including the SPHERE-IFS YH-band spec- 
trum (Cheetham et al. 2019), and SPHERE-IRDIS HB, 
K1, and K2-band photometry (Chauvin et al. 2017; 
Cheetham et al. 2019; also see Appendix D). 

To explore how the luminosity estimate depends on 
the choice of model spectrum used, we performed the 
analysis with a broad range of models spanning Tyg from 
1200 K to 1900 K and log(g) spanning 3.5 to 5.5. These 
models were drawn from three different grids, including 
two with different cloud implementations — BT-SETTL 
(Baraffe et al. 2015), and DRIFT-Phoenix (Witte et al. 
2009) — and the cloud-free Sonora-bobcat models (Mar- 
ley et al. 2021). No matter which model we used to fill 
in the gaps between the measured portions of the SED, 
log(Lpoi/Le) is always between —4.21 and —4.35. Con- 
sequently, the luminosity is constrained at the 0.15 dex 
level and the result is robust across all considered model 
atmospheres. 


4.6. Estimates of Mass and other Companion 
Properties from Hot-Start Evolutionary Models 


The mass of HIP 65426 b is estimated using a method 
similar to that described in Dupuy & Liu (2017). We 
first built an interpolated grid of model luminosities as a 
function of age and mass, with 10,000 equally-spaced age 
values spanning from 5 to 30 Myr and 10,000 equally- 
spaced mass values spanning from 0.3 to 21 Mjup using 
numpy.griddata in python . We adopted an age for 
HIP 65426 b of 1444 Myr based on the Lower Centaurus- 
Crux age given in Chauvin et al. (2017) and a mea- 
sured bolometric luminosity uniformly distributed be- 
tween log(Lpoi/Le) = —4.21 and —4.35 from the previ- 
ous section. 

We then generated 1x10° samples of age and mass 
from a Gaussian distribution in age around 14 Myr, with 
o=4 Myr, and a uniform distribution in mass from 0.3 
to 21 Myup. For each age, mass sample, we then look 
up the corresponding model luminosity from the inter- 
polated grid of model luminosities. For each age, mass 
sample, we accept the sample if the corresponding model 
luminosity is within the measured range of uniformly- 
distributed bolometric luminosities and reject the sam- 
ple if it lies outside this range. 

We implemented this procedure using the hybrid cloud 
grid from Saumon & Marley (2008). Given that this is 
a dusty, young red object, we expect the Saumon & 
Marley (2008) models, which take clouds into account, 
to provide the most reliable estimates of the properties 
of these objects among the model choices available. To 
sample the corresponding effective temperatures, sur- 
face gravities, and radii corresponding to our accepted 
mass values, we built interpolated grids of model effec- 
tive temperatures, surface gravities, and radii with the 
same spacing in age and mass as for the interpolated grid 
of luminosities, then looked up the corresponding prop- 
erty in the appropriate table for each accepted age, mass 
sample. Histograms of the final set of accepted masses, 
effective temperatures, surface gravities, and radii for 
each model are shown in Fig. 10. The best value of each 
parameter was taken as the median of the accepted dis- 
tribution, with error bars given by the 68% confidence 
interval as calculated from the histogram of each dis- 
tribution. We find a mass of 7.1+1.1 Mjup, radius of 
1.45+0.03 Ryup (which therefore imply a surface grav- 
ity log(g) = 3.937903), and effective temperature Teg = 
12822) Ki: 


4.7. Atmospheric Forward Model Comparison 


To explore the atmospheric properties of HIP 65426 b 
we performed a forward modelling analysis using the 
tool ForMoSA (Petrus et al. 2020), which compares spec- 
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Figure 10. Histograms of the final sets of accepted model properties for the hybrid cloud grid from Saumon & Marley (2008). 
The median value for each property is shown as a solid black line, with the 68% confidence region falling between the two dotted 


black lines. 


troscopic and/or photometric data with grids of precom- 
puted synthetic spectra. The code is based on the nested 
sampling algorithm (Skilling 2004), a Bayesian inversion 
method, that allows a global exploration of the param- 
eter space provided by the grid. In this work, we lim- 
ited our analysis to the BT-SETTL grid (CIFIST version, 
Allard et al. 2012; Baraffe et al. 2015) that accounts 
for convection using mixing-length theory, and works at 
hydrostatic, radiative-convective, and chemical equilib- 
rium. 

For the fit, we used a data set composed of the low 
resolution spectra (R, ~54) between 1.00 and 1.65 wm 
provided by VET /SPHERE-IFS (Chauvin et al. 2017), 
VLT /SPHERE-IRDIS Ap» (1.58 wm), H3 (1.66 wm), ky 
(2.11 um), and K2 (2.25 wm) photometry (Cheetham 
et al. 2019), and our new JWST NIRCam and MIRI 
photometry. The NaCo data are not used for the reasons 
described in Section 4.4. We first adapted the BT-SETTL 
synthetic spectra to our data, by reducing their spec- 
tral resolution to that of SPHERE-IFS and calculating 
the synthetic photometric flux at each bandpass using 
throughputs as obtained from spaceKLIP for the JWST 
data, and from the SVO filter service for all other data‘? 
(Rodrigo et al. 2012; Rodrigo & Solano 2020). We then 
defined flat priors on the Teg and the log(g) according 
to the limits of the grid, and applied nested sampling 
to estimate the posterior distributions of these two pa- 
rameters. We also add the radius, R, to the list of the 
parameters explored by the nested sampling. At each 
iteration, a radius is picked randomly (uniform prior), 
and a dilution factor Cx = (R/d)? is calculated and 
multiplied to the model, where d is the distance to the 
object (107.49 pc). 

The best fit models to our data combined with 
the existing SPHERE data are displayed in Fig- 


15 http://svo2.cab.inta-csic.es/theory /fps/ 


ure 9, alongside posteriors in Figure 11. We es- 
timate Teg=1667153 K, log(g)=4.0710 32 dex, and 
R=0.9210'}4 Ryup- From the Tyg and the radius, we ap- 
ply the Stefan-Boltzmann law and estimate a bolometric 
luminosity of log(Lpoi/Lo)=—4-23193, and from the 
log(g) and R we estimate a mass of M=3.95*7:78 Myup- 
By comparing the integrated flux of the best fit model 
across all wavelengths, and the integrated flux be- 
tween the shortest and longest wavelength observations, 
we determine that these observations span ~97% of 
HIP 65426 b’s luminous range. These results are also in 
agreement with a similar BT-SETTL model fitting proce- 
dure to VET /SINFONI data of HIP 65426 b (see Petrus 
et al. 2021). Finally, the uncertainties of all parameters 
are given at 30, however, we emphasise that they do not 
necessarily describe our confidence in the true planetary 
properties and are better considered as the model phase 
space that best fits our data. 

The precision on these measurements is primarily 
driven by the SPHERE/IFS data, however, we do 
note some differences that result from the addition of 
the JWST data. Specifically, when fitting just the 
SPHERE data in isolation, we obtain Teg = 1642733 K, 
log(g) = 3.89t3:1% dex, R = 1.0203) Ryup and 
log(Lpoi/Lo) = —4.10°30%, again with uncertainties 
given at 30. Therefore, the JWST data improve the 
precision of the radius and bolometric luminosity by a 
factor of ~3, but do not significantly improve the pre- 
cision on the temperature and surface gravity. We also 
note offsets between the best fit values between these two 
fits of up to 30 (using the uncertainties from the com- 
bined fit), which suggests that the selected models do 
not accurately capture these two regions of wavelength 
space simultaneously. 

The atmospheric forward model fit yields a luminos- 
ity within the bolometric luminosity range of —4.21 
to —4.35 that was found from combining SED mea- 
surements with models in the regions not covered by 
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Figure 11. Posterior distributions for the BT—Settl atmospheric model fitting to both JWST and VLT/SPHERE observations 
of HIP 65426b. Best fit values and lo uncertainties are indicated, however, these should be interpreted as the model phase 
space that fits these data, and not the precision to which these properties can be empirically measured. 


the SED. However, the effective temperatures and radii 
found using the atmospheric forward model are con- 
siderably in tension with predictions from evolution- 
ary model fits to the measured bolometric luminosity 
range (see Section 4.6). In particular, to obtain similar 
bolometric luminosity values, the forward models favour 
higher effective temperatures (~1700 K) and smaller 
radii (~0.90 Ryup) compared to the evolutionary mod- 
els (~1300 K, ~1.48 Ryup). In fact, the atmospheric 
forward model used here predicts an unphysically small, 
sub-Jupiter radius for an exoplanet which may still be 
contracting at these young ages. Thus, we consider the 
effective temperature and radius predictions from the 
evolutionary models to be more robust here. The ten- 
sion we find between the atmospheric models and evo- 
lutionary models is well-documented in the community; 
atmospheric models have a long history of requiring un- 


physically small radii and high effective temperatures to 
fit spectroscopy (see e.g., Marois et al. 2008; Patience 
et al. 2012). 


4.8. Future Work 


There is a range of additional investigation that can 
be performed on the data presented here that is worth 
highlighting, but ultimately falls outside the scope of 
this work. 

Most importantly, it is possible to perform the atmo- 
spheric forward model fitting procedure shown here in 
Section 4.7 with a wide range of state-of-the-art mod- 
els (e.g., ATMO, Tremblin et al. (2015); Phillips et al. 
(2020); Exo-REM, Charnay et al. (2018); Sonora, Kar- 
alidi et al. (2021)), each with their own treatment for 
the effects of clouds and atmospheric chemistry. Addi- 
tionally, atmospheric fitting can also be performed using 
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retrieval techniques (e.g., Molliére et al. 2020; Gonza- 
les et al. 2021). Divergences in the measured planetary 
properties between these models are expected, and a 
more complete analysis in the context of the relative 
assumptions of each model will greatly improve our un- 
derstanding of the true properties of HIP 65426 b. 

The precision of the 3—5yum data may be sufficient to 
provide constraints on the relative atmospheric abun- 
dances of CH, and CO, which can be impacted by dise- 
quilibrium chemistry (Zahnle & Marley 2014; Miles et al. 
2020). The F1140C photometry falls slightly under the 
best fit model at ~1.50, and may be indicative of ab- 
sorption by small silicate dust grains (Cushing et al. 
2006; Suarez & Metchev 2022, Miles et al. submitted). 
Finally, while the F1550C photometry is relatively in- 
sensitive to the choice of atmospheric model, it may be 
sensitive to circumplanetary disk emission (Sterzik et al. 
2004; Stolker et al. 2020a). 


5. CONCLUSION 


In this work we present the first ever scientific obser- 
vations using the JWST high contrast imaging modes of 
both NIRCam from 2—5 wm, and MIRI from 11—16 um. 
The known exoplanet companion, HIP 65426 b, is clearly 
detected in all seven observational filters, representing 
the first ever direct detection of an exoplanet beyond 
5 pm. These observations provide a variety of insights 
into: a) the performance and best practices of JWST 
high contrast imaging, and b) the properties of the 
HIP 65426 b system, which we summarise below: 


e Contrast: JWST is exceeding its anticipated 
contrast performance for both NIRCam and MIRI 
coronagraphy by up to a factor of 10 in the con- 
trast limited regime (see Section 4.1). For the 
contrasts achieved, we are sensitive to sub-Jupiter 
companions with masses as small as 0.3Mjup 
beyond separations of ~100 au. Furthermore, 
for more optimal targets such as young, nearby 
M stars it is highly likely that both NIRCam and 
MIRI will be sensitive to sub-Saturn mass objects 
beyond ~10 au (Carter et al. 2021a). 


e Subtraction Strategy: For these data at small 
separations <2”, the best contrast is obtained 
using a small-grid dither RDI subtraction strat- 
egy for both NIRCam and MIRI. Additionally, an 
ADI+RDI subtraction does not significantly im- 
prove the measured contrast compared to RDI. 
For the MIRI F1550C observations in particular, 
we were unable to recover HIP 65426 b using ADI 
alone. At wider separations however, the observa- 
tional efficiency of ADI may make it preferable to 


RDI. These conclusions may aid future observers 
in selecting their PSF subtraction strategy, al- 
though we emphasise that a clearer understanding 
of whether they apply under all circumstances will 
require the analysis of a broader range of JWST 
coronagraphic observations. 


Photometry: These photometric observations of 
HIP 65426 b provide exquisite sensitivity at a pre- 
cision of ~7% for NIRCam and ~16% for MIRI. 
Furthermore, prior to propagation of the uncer- 
tainty in the stellar flux (~5%), and the current 
absolute flux calibration accuracy (5/15% for NIR- 
Cam/MIRI), the uncertainty in the measured rela- 
tive flux is even smaller at ~2% for both NIRCam 
and MIRI. These measurements are a significant 
step forward from ground-based observations from 
3—5 pm, which have comparative sensitivities of 
~13—32% for HIP 65426 b, and are restricted to 
particular wavelength regions due to telluric con- 
tamination. With this improved precision we will 
be able to constrain directly imaged exoplanet at- 
mospheres in much greater detail, in addition to 
more complex effects such as variability, disequi- 
librium chemistry, and the emission of circumplan- 
etary material. 


Atmospheric Model Fitting: Using a 
BT-SETTL atmospheric forward model we are 
able to fit all data, in addition to the 
majority of ground-based observations to 
within 20. This agreement provides pre- 
cise constraints on the Tog = 1667*3?) K, 
log(g) = 4.077513 dex, R = 0.927004 Ryup, and 
log(Lyoi/Lo)=—4.231)'35. Compared to a fit ex- 
cluding the JWST data, this corresponds to a 
factor of ~3 improvement in the precision of the 
radius and bolometric luminosity. Despite the 
excellent model agreement, both the temperature 
and unphysically small radius are in disagreement 
with the values obtained from the evolutionary 
models, further emphasising a long standing ten- 
sion for this class of objects. 


Empirical Bolometric Luminosity: As JWST 
offers a uniquely broad spectral coverage in com- 
parison to ground-based instruments, we are able 
to obtain a very precise measurement of the bolo- 
metric luminosity of HIP65426b that is con- 
strained between a log(Lyo1/Lo)=—4.21 to —4.35, 
irrespective of the model atmosphere adopted 
for the wavelengths not covered by observa- 
tions. In combination with evolutionary mod- 
els, this provides tight constraints on the prop- 
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erties of HIP65426b with M=7.141.1 Myup, 
Topp =128213° K, and R=1.45+0.03 Ryup. Given 
the achieved sensitivity, similar JWST observa- 
tions will facilitate this analysis for a broader 
range of PMCs than ever before and provide com- 
parable constraints on their bolometric luminosi- 
ties, and therefore mass. 


In conclusion, the observations reported here from our 
ERS program (ERS-01386, Hinkley et al. 2022a) demon- 
strate that JWST provides an transformative opportu- 
nity to study exoplanets through high contrast imaging. 
Beyond this work, we look forward to further results 
from our program of: 3-16 wm NIRCam and MIRI coro- 
nagraphy of the circumstellar disk HD 141569 A (Millar- 
Blanchaer et al. in preparation, Choquet et al. in prepa- 
ration); NIRSpec and MIRI spectroscopy from 1—28 zm 
of the PMC VHS J1256b (Miles et al. submitted); and 
NIRISS AMI observations of HIP 65426 at 3.8 wm (Ray 
et al. in preparation, Sallum et al. in preparation). 
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APPENDIX 
A. SUBTRACTED IMAGES 


(4s/Afin) S}UNOD 


Figure 12. As in Figure 3, except for the NIRCam F250M, F300M, F356W, F410M, and F444W filters. Here we show 


subtractions using the maximum number of PCA modes for ADI, RDI, and ADI+RDI, respectively. The F250M subtracted 


images have been smoothed as described in Section 2. 
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Figure 13. As in Figure 3, except for the MIRI F1140C and F1550C filters. Here we show subtractions using the maximum 
number of PCA modes for ADI, RDI, and ADI+RDI, respectively. To aid visual clarity, the subtracted F1550C images are 
shown with a peak image intensity five times smaller than the unsubtracted image. 
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B. CONTRAST PERFORMANCE 
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Figure 14. As in Figure 5, but for all used filters. We also plot the equivalent predicted contrast curves for these observations 
from PanCAKE (grey lines) following (Carter et al. 2021b). In every filter JWST is exceeding its predicted performance. 
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C. PSF FITTING 
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Figure 15. The data (left column), model PSF (middle column), and residuals (right column) for the spaceKLIP PSF fitting 
of the NIRCam observations of HIP 65426 b. Pixel counts are in MJy/sr and are indicated by the colour bar on the right hand 


side. 
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Figure 16. As in Figure 15, but for the MIRI observations. 
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D. COMPLEMENTARY PHOTOMETRIC MEASUREMENTS OF HIP 65426 B 


Table 3. Additional photometric measurements of HIP 65426 b considered in this work. 


Wavelength (um) Bandwidth (wm) Instrument Band Flux (W m~?ym7?) Ref. 
1.002 0.011 SPHERE/IFS YJ (2.434£0.569)x 10-271 
1.011 0.011 SPHERE/IFS YJ (3.155 + 0.670) x 107271 
1.021 0.011 SPHERE/IFS YJ (3.564 0.575) x 107271 
1.030 0.011 SPHERE/IFS YJ (3.199 0.428) x 107271 
1.040 0.011 SPHERE/IFS YJ (4.095 + 0.482) x 107271 
1.050 0.011 SPHERE/IFS YJ (3.829 0.435) x 107271 
1.060 0.011 SPHERE/IFS YJ (4.154£0.489)x 107171 
1.070 0.011 SPHERE/IFS YJ (4.456-£0.521)x 107171 
1.081 0.011 SPHERE/IFS YJ (4.679 0.432) x 107171 
1.091 0.011 SPHERE/IFS YJ (5.1434 0.483) x 107271 
1.102 0.011 SPHERE/IFS YJ (4.746£0.543)x 107271 
1.112 0.011 SPHERE/IFS YJ (5.1274 0.576) x 10-71 
1.123 0.011 SPHERE/IFS YJ (4.775 40.641) x 107271 
1.133 0.011 SPHERE/IFS YJ (5.2284 0.682) x 107271 
1.144 0.011 SPHERE/IFS YJ (4.9824 0.609) x 10771 
1.154 0.011 SPHERE/IFS YJ (4.82140.519)x 107271 
1.165 0.011 SPHERE/IFS YJ (5.32440.486)x 107271 
1.175 0.011 SPHERE/IFS YJ (4.7394 0.453) x 10771 
1.186 0.011 SPHERE/IFS YJ (5.7094 0.477) x 10771 
1.196 0.011 SPHERE/IFS YJ (5.29140.415)x 107271 
1.206 0.011 SPHERE/IFS YJ (6.15540.471)x 107271 
1.217 0.011 SPHERE/IFS YJ (6.586+0.490) x 10771. 
1.227 0.011 SPHERE/IFS YJ (6.43140.500) x 10-71 
1.237 0.011 SPHERE/IFS YJ (6.2034 0.483) x 10771 
1.247 0.011 SPHERE/IFS YJ (6.40940.477)x 10-71 
1.257 0.011 SPHERE/IFS YJ (7.1004 0.483) x 10771 
1.266 0.011 SPHERE/IFS YJ (7.289+0.502)x 10771 
1.276 0.011 SPHERE/IFS YJ (7.344£0.504)x 107271 
1.285 0.011 SPHERE/IFS YJ (7.31540.499)x 10-71 
1.294 0.011 SPHERE/IFS YJ (8.21540.543)x 107271 
1.303 0.011 SPHERE/IFS YJ (8.70140.543)x 107271 
1.312 0.011 SPHERE/IFS YJ (8.731£0.580) x 107271 
1.321 0.011 SPHERE/IFS YJ (7.45040.594) x 10771 
1.329 0.011 SPHERE/IFS YJ (7.073 0.583) x 107271 
0.987 0.019 SPHERE/IFS YJH  (1.55640.312)x 10-17 1 
1.002 0.019 SPHERE/IFS YJH = (1.79140.355)x107!7 1 
1.018 0.019 SPHERE/IFS YJH  (2.855+0.510)x 107171 
1.034 0.019 SPHERE/IFS YJH  (2.615+0.305)x 10717 1 
1.051 0.019 SPHERE/IFS YJH  (3.394+0.442)x 107171 
1.068 0.019 SPHERE/IFS YJH  (4.239+0.456)x 107271 
1.086 0.019 SPHERE/IFS YJH (3.472+0.407)x 107271 
1.104 0.019 SPHERE/IFS YJH (3.802+0.415)x 107171 
1.122 0.019 SPHERE/IFS YJH (4.403+0.504)x 10-171 


Table 3 continued 


JWST IMAGING OF HIP 65426B FROM 2—16 MICRON 
Table 3 (continued) 


Wavelength (um) Bandwidth (wm) Instrument Band Flux (W m~?ym7?) Ref. 
1.140 0.019 SPHERE/IFS YJH  (3.906+0.397)x 107!" 1 
1.159 0.019 SPHERE/IFS YJH  (4.18240.422)x 1077 1 
1.178 0.019 SPHERE/IFS YJH  (4.64140.452)x 1077 1 
1.197 0.019 SPHERE/IFS YJH = (5.868 +0.544)x 107'7 1 
1.216 0.019 SPHERE/IFS YJH  (6.356+0.573)x 107!" 1 
1.235 0.019 SPHERE/IFS YJH = (6.7634 0.602) x 107'7 1 
1.255 0.019 SPHERE/IFS YJH = (7.107 +0.622)x 10-71 
1.274 0.019 SPHERE/IFS YJH = (7.228+0.630)x 10717 1 
1.294 0.019 SPHERE/IFS YJH (7.599 + 0.668) x 10°!" 1 
1.314 0.019 SPHERE/IFS YJH = (7.296 + 0.648) x 1077 1 
1.333 0.019 SPHERE/IFS YJH = (6.046+40.577)x 107'7 1 
1.353 0.019 SPHERE/IFS YJH  (5.400+0.649)x 1077 1 
1.372 0.019 SPHERE/IFS YJH = (5.523+0.872)x 107" 1 
1.391 0.019 SPHERE/IFS YJH  (5.875+0.729)x 107!" 1 
1.411 0.019 SPHERE/IFS YJH (5.285 +0.550)x 10771 
1.430 0.019 SPHERE/IFS YJH  (4.516+0.460)x 107’ 1 
1.449 0.019 SPHERE/IFS YJH (4.768 + 0.439) x 107'7 1 
1.467 0.019 SPHERE/IFS YJH = (5.874+0.484)x 1077 11 
1.486 0.019 SPHERE/IFS YJH (5.888 +0.506)x 107'7 1 
1.504 0.019 SPHERE/IFS YJH  (5.92740.498)x 1077 1 
1.522 0.019 SPHERE/IFS YJH = (6.259+ 0.519) x 107" 1 
1.539 0.019 SPHERE/IFS YJH = (6.778 40.561) x 1077 1 
1.556 0.019 SPHERE/IFS YJH = (7.318 40.605) x 107'7 1 
1.573 0.019 SPHERE/IFS YJH = (7.560 + 0.624) x 107" 1 
1.589 0.019 SPHERE/IFS YJH (8.037 + 0.665) x 107" 1 
1.605 0.019 SPHERE/IFS YJH  (8.497+0.706) x 107" 1 
1.621 0.019 SPHERE/IFS YJH = (8.579+0.715)x 1077 1 
1.636 0.019 SPHERE/IFS YJH  (9.01140.762)x 1077 1 
1.593 0.055 SPHERE/IRDIS  H2 (8.569 + 0.383) x 10717 1 
1.667 0.056 SPHERE/IRDIS  4H3 (10.129+ 0.564) x 107*7 1 
2.110 0.102 SPHERE/IRDIS = K1 (7.500 + 0.600) x 10717 1 
2.251 0.109 SPHERE/IRDIS  K2 (7.100 + 0.600) x 1077 1 
3.800 0.620 NACO L’ (4.010 + 0.542) x 1071-2, 3 
4.051 0.020 NACO NB4.05 (3.220+0.780) x 10717 3 
4.780 0.590 NACO M’ (2.549 + 0.820) x 10717) -2,3 


References—(1) Chauvin et al. (2017); (2) Cheetham et al. (2019); (3) Stolker et al. (2020b). 


Note—For the L’ and M’-band photometry, we considered the average of the measurements reported in 
Cheetham et al. (2019) and Stolker et al. (2020b), but kept the largest error bars. 


